Notebook

R Markdown

This is an R Markdown document. Markdown is a simple formatting syntax for authoring HTML, PDF, and MS Word documents. For more details on using R Markdown see http://rmarkdown.rstudio.com.

When you click the Knit button a document will be generated that includes both content as well as the output of any embedded R code chunks within the document. You can embed an R code chunk like this:

install.packages("dplyr")

Error in install.packages : Updating loaded packages

library(dplyr)

## report mean ki per ROI

dataMNI <- read.csv('~/Dropbox/DCCN/Creativity/roi_results_mni.csv')
  dataMNI %>%
  group_by(ROI) %>%
  dplyr::summarise(mean = mean(meanRoiContrastEstimate),
  sd = sd(meanRoiContrastEstimate)) 

NA

install.packages("tidyverse")

Error in install.packages : Updating loaded packages

library(tidyverse)

## graph the distribution in box plot

ggplot(dataMNI, aes(x=ROI, y=meanRoiContrastEstimate)) + 
    geom_boxplot() + 
    coord_flip()

## look at the density plot

ggplot(dataMNI, aes(x=meanRoiContrastEstimate, group=ROI)) + 
    geom_line(stat="density")

## report mean ki per ROI

dataNative <- read.csv('~/Dropbox/DCCN/Creativity/roi_results_native.csv')
install.packages("dplyr")

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.0/dplyr_1.0.4.tgz'
Content type 'application/x-gzip' length 1250919 bytes (1.2 MB)
==================================================
downloaded 1.2 MB

The downloaded binary packages are in
    /var/folders/wq/bqdw5b953sbc3lfdvxcslm2c0000gn/T//RtmpMcvwGX/downloaded_packages

install.packages("tidyverse")

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.0/tidyverse_1.3.0.tgz'
Content type 'application/x-gzip' length 433049 bytes (422 KB)
==================================================
downloaded 422 KB

The downloaded binary packages are in
    /var/folders/wq/bqdw5b953sbc3lfdvxcslm2c0000gn/T//RtmpMcvwGX/downloaded_packages

  dataNative %>%
  group_by(ROI) %>%
  dplyr::summarise(mean = mean(meanRoiContrastEstimate),
  sd = sd(meanRoiContrastEstimate)) 

NA

## graph the distribution in box plot

ggplot(dataNative, aes(x=ROI, y=meanRoiContrastEstimate)) + 
    geom_boxplot() + 
    coord_flip()

## look at the density plot

ggplot(dataNative, aes(x=meanRoiContrastEstimate, group=ROI)) + 
    geom_line(stat="density")

## Ruben shared that we are only using native space and only whole caudate, putamen and VS
install.packages("sjPlot")

Error in install.packages : Updating loaded packages

library(sjPlot)
dataAll <- read.csv('~/Dropbox/DCCN/Creativity/CreativityIDsheet_MissingDataRemoved.csv')

# turn drug conditions into factor levels
dataAll$Drug <- factor(dataAll$Drug, levels = c("MPH","SUL","PBO"));

# set contrasts to sum-to-zero
options(contrasts=c("contr.sum", "contr.poly"))

# set two dataframes for the contrast between MPH and PBO - between SUL and PBO
df_MPH <- dplyr::filter(dataAll, Drug %in% c("MPH","PBO"))  
df_SUL <- dplyr::filter(dataAll, Drug %in% c("SUL","PBO"))

SelectData <- dataNative[grepl("wholeCaudate", dataNative[["ROI"]]) | grepl("wholePutamen", dataNative[["ROI"]]) | grepl("VS", dataNative[["ROI"]]), ]

SelectData <- SelectData[-grep("_", SelectData$ROI),]

ggplot(SelectData, aes(x=ROI, y=meanRoiContrastEstimate)) + 
    geom_boxplot() + 
    coord_flip()

ggplot(SelectData, aes(x=meanRoiContrastEstimate, group=ROI)) + 
    geom_line(stat="density")

# How do different ROI Ki s relate to each other?
my_cols <- c("#00AFBB", "#E7B800", "#FC4E07")  
pairs(dataAll[,2:4], pch = 19,  cex = 0.5,
      lower.panel=NULL)

# Does session number have an effect on convergent thinking?

ggplot(dataAll, aes(x=as.factor(Session), y=Convergent_Pasta)) + 
    geom_boxplot(fill="slateblue", alpha=0.2) + 
    xlab("Session")

# Total number of generated ideas

SessionEffectTotal <- lmer(Total ~ 1 + Session + (1 + Session|| ID), data = dataAll, REML=F) 
print(summary(SessionEffectTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Total ~ 1 + Session + ((1 | ID) + (0 + Session | ID))
   Data: dataAll

     AIC      BIC   logLik deviance df.resid 
  2202.0   2220.2  -1096.0   2192.0      277 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-5.0605 -0.4235 -0.0551  0.3522  4.1868 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 167.328  12.936  
 ID.1     Session       9.404   3.067  
 Residual              56.000   7.483  
Number of obs: 282, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)  28.0887     1.7805  15.776
Session       2.3457     0.6308   3.719

plot_model(SessionEffectTotal, type = "pred", show.data=TRUE, terms = c("Session"))
install.packages("sjPlot")

trying URL 'https://cran.rstudio.com/bin/macosx/contrib/4.0/sjPlot_2.8.7.tgz'
Content type 'application/x-gzip' length 1499146 bytes (1.4 MB)
==================================================
downloaded 1.4 MB

The downloaded binary packages are in
    /var/folders/wq/bqdw5b953sbc3lfdvxcslm2c0000gn/T//RtmpMcvwGX/downloaded_packages

# Pure convergent scores

SessionEffect <- lmer(Convergent_Pasta ~ 1 + Session + (1 + Session|| ID), data = dataAll, REML=F) 
print(summary(SessionEffect), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Convergent_Pasta ~ 1 + Session + ((1 | ID) + (0 + Session | ID))
   Data: dataAll

     AIC      BIC   logLik deviance df.resid 
  2183.2   2201.4  -1086.6   2173.2      277 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-4.1205 -0.4193 -0.0678  0.3707  5.0262 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 144.02   12.001  
 ID.1     Session       9.27    3.045  
 Residual              53.67    7.326  
Number of obs: 282, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)  16.9078     1.6924   9.990
Session       2.8777     0.6198   4.643

plot_model(SessionEffect, type = "pred", show.data=TRUE, terms = c("Session"))

# convergent/divergent ratio

SessionEffectRatio <- lmer(Con_Div ~ 1 + Session + (1 + Session|| ID), data = dataAll, REML=F) 
print(summary(SessionEffectRatio), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Div ~ 1 + Session + ((1 | ID) + (0 + Session | ID))
   Data: dataAll

     AIC      BIC   logLik deviance df.resid 
  1796.2   1814.1   -893.1   1786.2      260 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.0127 -0.2844 -0.2047 -0.0056  5.2673 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept)  9.917   3.149   
 ID.1     Session      3.854   1.963   
 Residual             31.422   5.606   
Number of obs: 265, groups:  ID, 93

Fixed effects:
            Estimate Std. Error t value
(Intercept)   3.0677     0.9668   3.173
Session       0.9149     0.4753   1.925

plot_model(SessionEffectRatio, type = "pred", show.data=TRUE,terms = c("Session"))

NA
NA

# Does session number have an effect on divergent thinking?

ggplot(dataAll, aes(x=as.factor(Session), y=Divergent_Pasta)) + 
    geom_boxplot(fill="slateblue", alpha=0.2) + 
    xlab("Session")

SessionEffect2 <- lmer(Divergent_Pasta ~ 1 + Session + (1 + Session|| ID), data = dataAll, REML=F) 
print(summary(SessionEffect2), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Divergent_Pasta ~ 1 + Session + ((1 | ID) + (0 + Session | ID))
   Data: dataAll

     AIC      BIC   logLik deviance df.resid 
  1872.3   1890.5   -931.1   1862.3      277 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.2533 -0.5594 -0.1156  0.4431  3.5352 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 27.342   5.229   
 ID.1     Session      2.999   1.732   
 Residual             21.882   4.678   
Number of obs: 282, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)  11.1809     0.9133  12.243
Session      -0.5319     0.3851  -1.381

plot_model(SessionEffect2, type = "pred", show.data=TRUE, terms = c("Session"))

## look at the correlation between convergent and divergent thinking for each session
my_cols <- c("#00AFBB", "#E7B800", "#FC4E07")  
pairs(dataAll[,7:8], pch = 19,  cex = 0.5,
      col = my_cols[dataAll$Session],
      lower.panel=NULL)

# Does drug on its own have an effect on Convergent thinking?

ggplot(dataAll, aes(x=as.factor(Drug), y=Convergent_Pasta)) + 
    geom_boxplot(fill="slateblue", alpha=0.2) + 
    xlab("Session")

DrugEffect <- lmer(Convergent_Pasta ~ 1 + Drug + (1 | ID), data = df_MPH, REML=F) 
print(summary(DrugEffect), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Convergent_Pasta ~ 1 + Drug + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1523.8   1536.8   -757.9   1515.8      184 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.8151 -0.4594 -0.0832  0.3814  3.7880 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 159.59   12.633  
 Residual              85.41    9.242  
Number of obs: 188, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)  22.7447     1.4670  15.504
Drug1         0.3085     0.6740   0.458

plot_model(DrugEffect, type = "pred", show.data=TRUE, terms = c("Drug"))

# Drug effect in the presence of session?

DrugAndSessionEffect <- lmer(Convergent_Pasta ~ 1 + Drug*Session + (1 | ID), data = df_MPH, REML=F) 
print(summary(DrugAndSessionEffect), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Convergent_Pasta ~ 1 + Drug * Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1520.3   1539.7   -754.1   1508.3      182 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.2378 -0.4158 -0.0723  0.3572  4.1477 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 165.74   12.874  
 Residual              77.86    8.824  
Number of obs: 188, groups:  ID, 94

Fixed effects:
              Estimate Std. Error t value
(Intercept)    18.1707     2.2385   8.117
Drug1          -1.8587     2.5543  -0.728
Session         2.3334     0.8697   2.683
Drug1:Session   1.0552     1.2772   0.826

plot_model(DrugAndSessionEffect, type = "pred", show.data=TRUE, terms = c("Drug","Session"))

# Does drug on its own have an effect on Divergent thinking?

ggplot(dataAll, aes(x=as.factor(Drug), y=Divergent_Pasta)) + 
    geom_boxplot(fill="slateblue", alpha=0.2) + 
    xlab("Session")

DrugEffect2 <- lmer(Divergent_Pasta ~ 1 + Drug + (1 | ID), data = df_MPH, REML=F) 
print(summary(DrugEffect2), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Divergent_Pasta ~ 1 + Drug + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1290.9   1303.8   -641.4   1282.9      184 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.3879 -0.4947 -0.0901  0.3593  3.6678 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 37.76    6.145   
 Residual             28.00    5.292   
Number of obs: 188, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)  10.1862     0.7421  13.727
Drug1         0.6968     0.3859   1.805

plot_model(DrugEffect2, type = "pred", terms = c("Drug"))

# Drug effect in the presence of session?

DrugAndSessionEffect2 <- lmer(Divergent_Pasta ~ 1 + Drug* Session + (1 | ID), data = df_MPH, REML=F) 
print(summary(DrugAndSessionEffect2), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Divergent_Pasta ~ 1 + Drug * Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1292.0   1311.4   -640.0   1280.0      182 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.2919 -0.4825 -0.0912  0.3654  3.7287 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 36.27    6.023   
 Residual             27.97    5.289   
Number of obs: 188, groups:  ID, 94

Fixed effects:
              Estimate Std. Error t value
(Intercept)    11.1113     1.2360   8.990
Drug1           2.6937     1.4175   1.900
Session        -0.4498     0.5148  -0.874
Drug1:Session  -1.0190     0.7047  -1.446

plot_model(DrugAndSessionEffect2, type = "pred", terms = c("Drug","Session"))

# How do different ROI Ki s relate to each other?
cols <- c(2:4, 7:8)
my_cols <- c("#00AFBB", "#E7B800", "#FC4E07")  
pairs(dataAll[,cols], pch = 19,  cex = 0.5,
            col = my_cols[dataAll$Session],
      lower.panel=NULL)

############################# Test the effects of Caudate Ki #############################################################

## effects of MPH vs PBO below

# pure Convergent score
MPH_Caudate_PureConvergent <- lmer(Convergent_Pasta ~ 1 + Drug*Caudate_ki + Session + (1 | ID), data = df_MPH, REML=F) 
print(summary(MPH_Caudate_PureConvergent), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Convergent_Pasta ~ 1 + Drug * Caudate_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1522.0   1544.7   -754.0   1508.0      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.1830 -0.4290 -0.0762  0.3575  4.1319 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 165.85   12.878  
 Residual              77.65    8.812  
Number of obs: 188, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)       18.7640    10.3977   1.805
Drug1             -4.0905     4.4653  -0.916
Caudate_ki       -50.2712   722.7969  -0.070
Session            2.4191     0.8713   2.776
Drug1:Caudate_ki 305.3485   315.6846   0.967

plot_model(MPH_Caudate_PureConvergent, type = "pred", terms = c("Caudate_ki","Session","Drug"))

plot(MPH_Caudate_PureConvergent)  # Plot the model information

qqnorm(residuals(MPH_Caudate_PureConvergent)) 

# Convergent/Divergent ratio score
MPH_Caudate_ConDivRatio <- lmer(Con_Div ~ 1 + Drug*Caudate_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_Caudate_ConDivRatio), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Div ~ 1 + Drug * Caudate_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1226.7   1248.8   -606.3   1212.7      168 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.7803 -0.3750 -0.2319 -0.0082  5.1290 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 20.16    4.490   
 Residual             42.87    6.547   
Number of obs: 175, groups:  ID, 91

Fixed effects:
                 Estimate Std. Error t value
(Intercept)        3.6777     4.8746   0.754
Drug1             -3.5583     3.3911  -1.049
Caudate_ki       -27.9163   330.1346  -0.085
Session            0.8616     0.6456   1.335
Drug1:Caudate_ki 253.1895   239.3929   1.058

plot_model(MPH_Caudate_ConDivRatio, type = "pred", terms = c("Caudate_ki","Session","Drug"))

plot(MPH_Caudate_ConDivRatio)

qqnorm(residuals(MPH_Caudate_ConDivRatio))

# Convergent/Total ratio score
MPH_Caudate_ConTotal <- lmer(Con_Total ~ 1 + Drug*Caudate_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_Caudate_ConTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Total ~ 1 + Drug * Caudate_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
   -54.5    -31.9     34.3    -68.5      179 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.0574 -0.5147  0.1010  0.5240  2.1860 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03478  0.1865  
 Residual             0.01852  0.1361  
Number of obs: 186, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)       0.50870    0.15288   3.327
Drug1            -0.05650    0.06897  -0.819
Caudate_ki        6.22253   10.60444   0.587
Session           0.03407    0.01357   2.511
Drug1:Caudate_ki  3.18761    4.87762   0.654

plot_model(MPH_Caudate_ConTotal, type = "pred", terms = c("Caudate_ki","Session","Drug"))

plot(MPH_Caudate_ConTotal)

qqnorm(residuals(MPH_Caudate_ConTotal))

# Convergent-Divergent difference score 
MPH_Caudate_ConDifference <- lmer(ConDiv_Dif ~ 1 + Drug*Caudate_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_Caudate_ConDifference), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: ConDiv_Dif ~ 1 + Drug * Caudate_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1609.7   1632.4   -797.9   1595.7      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.1195 -0.4122 -0.0729  0.4347  3.9454 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 232.5    15.25   
 Residual             134.8    11.61   
Number of obs: 188, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)         6.365     12.624   0.504
Drug1              -3.333      5.883  -0.566
Caudate_ki         56.472    874.925   0.065
Session             2.791      1.142   2.443
Drug1:Caudate_ki  199.909    415.911   0.481

plot_model(MPH_Caudate_ConDifference, type = "pred", terms = c("Caudate_ki","Session","Drug"))

plot(MPH_Caudate_ConDifference)

qqnorm(residuals(MPH_Caudate_ConDifference))

# Divergent/Total ratio score
MPH_Caudate_DivTotal <- lmer(Div_Total ~ 1 + Drug*Caudate_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_Caudate_DivTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Div_Total ~ 1 + Drug * Caudate_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
   -54.4    -31.8     34.2    -68.4      179 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.1847 -0.5250 -0.1048  0.5143  2.0560 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03477  0.1865  
 Residual             0.01855  0.1362  
Number of obs: 186, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)       0.49187    0.15288   3.217
Drug1             0.05697    0.06903   0.825
Caudate_ki       -6.22980   10.60419  -0.587
Session          -0.03418    0.01358  -2.517
Drug1:Caudate_ki -3.21665    4.88194  -0.659

plot_model(MPH_Caudate_DivTotal, type = "pred", terms = c("Caudate_ki","Session","Drug"))

plot(MPH_Caudate_DivTotal)

qqnorm(residuals(MPH_Caudate_DivTotal))

## effects of SUL vs PBO below

# pure Convergent score
SUL_Caudate_PureConvergent <- lmer(Convergent_Pasta ~ 1 + Drug*Caudate_ki + Session + (1 | ID), data = df_SUL, REML=F) 
print(summary(SUL_Caudate_PureConvergent), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Convergent_Pasta ~ 1 + Drug * Caudate_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
  1494.8   1517.5   -740.4   1480.8      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.9468 -0.4408 -0.0731  0.2915  3.9815 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 178.63   13.365  
 Residual              57.41    7.577  
Number of obs: 188, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)       18.1860    10.3607   1.755
Drug1             -4.0681     3.8259  -1.063
Caudate_ki       -88.7444   727.4964  -0.122
Session            2.7543     0.7906   3.484
Drug1:Caudate_ki 269.0754   270.9883   0.993

plot_model(SUL_Caudate_PureConvergent, type = "pred", terms = c("Caudate_ki","Session","Drug"))

plot(SUL_Caudate_PureConvergent)

qqnorm(residuals(SUL_Caudate_PureConvergent))

# Convergent/Divergent ratio score
SUL_Caudate_ConDivRatio <- lmer(Con_Div ~ 1 + Drug*Caudate_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_Caudate_ConDivRatio), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Div ~ 1 + Drug * Caudate_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
  1206.4   1228.7   -596.2   1192.4      171 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.8341 -0.3568 -0.2152  0.0065  5.7767 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 16.10    4.012   
 Residual             33.96    5.827   
Number of obs: 178, groups:  ID, 93

Fixed effects:
                  Estimate Std. Error t value
(Intercept)         8.3899     4.2775   1.961
Drug1               0.4849     3.0083   0.161
Caudate_ki       -369.9599   294.5373  -1.256
Session             0.6971     0.5930   1.175
Drug1:Caudate_ki  -61.0075   213.5831  -0.286

plot_model(SUL_Caudate_ConDivRatio, type = "pred", terms = c("Caudate_ki","Session","Drug"))

plot(SUL_Caudate_ConDivRatio)

qqnorm(residuals(SUL_Caudate_ConDivRatio))

# Convergent/Total ratio score
SUL_Caudate_ConTotal <- lmer(Con_Total ~ 1 + Drug*Caudate_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_Caudate_ConTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Total ~ 1 + Drug * Caudate_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
   -50.9    -28.2     32.4    -64.9      181 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.22401 -0.55586  0.07603  0.51605  1.83188 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03359  0.1833  
 Residual             0.01977  0.1406  
Number of obs: 188, groups:  ID, 94

Fixed effects:
                  Estimate Std. Error t value
(Intercept)       0.497809   0.150729   3.303
Drug1            -0.006216   0.071004  -0.088
Caudate_ki        1.950943  10.536689   0.185
Session           0.068229   0.014451   4.721
Drug1:Caudate_ki -0.859790   5.029105  -0.171

plot_model(SUL_Caudate_ConTotal, type = "pred", terms = c("Caudate_ki","Session","Drug"))

plot(SUL_Caudate_ConTotal)

qqnorm(residuals(SUL_Caudate_ConTotal))

# Convergent-Divergent difference score 
SUL_Caudate_ConDifference <- lmer(ConDiv_Dif ~ 1 + Drug*Caudate_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_Caudate_ConDifference), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: ConDiv_Dif ~ 1 + Drug * Caudate_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
  1589.6   1612.3   -787.8   1575.6      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.7010 -0.4469 -0.0126  0.3783  3.7805 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 217.3    14.74   
 Residual             118.1    10.87   
Number of obs: 188, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)         5.339     12.009   0.445
Drug1              -3.695      5.487  -0.673
Caudate_ki         75.444    840.056   0.090
Session             3.162      1.119   2.825
Drug1:Caudate_ki  221.314    388.615   0.569

plot_model(SUL_Caudate_ConDifference, type = "pred", terms = c("Caudate_ki","Session","Drug"))

plot(SUL_Caudate_ConDifference)

qqnorm(residuals(SUL_Caudate_ConDifference))

# Divergent/Total ratio score
SUL_Caudate_DivTotal <- lmer(Div_Total ~ 1 + Drug*Caudate_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_Caudate_DivTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Div_Total ~ 1 + Drug * Caudate_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
   -50.9    -28.3     32.5    -64.9      181 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-1.83278 -0.51699 -0.07626  0.55621  2.22345 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03356  0.1832  
 Residual             0.01977  0.1406  
Number of obs: 188, groups:  ID, 94

Fixed effects:
                  Estimate Std. Error t value
(Intercept)       0.502142   0.150691   3.332
Drug1             0.006258   0.071005   0.088
Caudate_ki       -1.940537  10.533992  -0.184
Session          -0.068224   0.014451  -4.721
Drug1:Caudate_ki  0.849165   5.029166   0.169

plot_model(SUL_Caudate_DivTotal, type = "pred", terms = c("Caudate_ki","Session","Drug"))

plot(SUL_Caudate_DivTotal)

qqnorm(residuals(SUL_Caudate_DivTotal))

NA
NA

############################# Test the effects of Putamen Ki #############################################################

## effects of MPH vs PBO below

# pure Convergent score
MPH_Putamen_PureConvergent <- lmer(Convergent_Pasta ~ 1 + Drug*Putamen_ki + Session + (1 | ID), data = df_MPH, REML=F) 
print(summary(MPH_Putamen_PureConvergent), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Convergent_Pasta ~ 1 + Drug * Putamen_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1522.9   1545.6   -754.5   1508.9      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.1949 -0.4529 -0.0786  0.3663  4.1514 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 165.35   12.859  
 Residual              78.46    8.858  
Number of obs: 188, groups:  ID, 94

Fixed effects:
                  Estimate Std. Error t value
(Intercept)       19.36174   11.43444   1.693
Drug1              0.04279    5.02649   0.009
Putamen_ki       -69.91372  664.58081  -0.105
Session            2.35409    0.88805   2.651
Drug1:Putamen_ki   8.36032  296.07382   0.028

plot_model(MPH_Putamen_PureConvergent, type = "pred", terms = c("Putamen_ki","Session","Drug"))

# Convergent/Divergent ratio score
MPH_Putamen_ConDivRatio <- lmer(Con_Div ~ 1 + Drug*Putamen_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_Putamen_ConDivRatio), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Div ~ 1 + Drug * Putamen_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1227.6   1249.8   -606.8   1213.6      168 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.8090 -0.3680 -0.2338 -0.0256  5.1362 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 20.20    4.495   
 Residual             43.15    6.569   
Number of obs: 175, groups:  ID, 91

Fixed effects:
                 Estimate Std. Error t value
(Intercept)        4.1489     5.3267   0.779
Drug1              1.4390     3.7759   0.381
Putamen_ki       -40.7150   303.5188  -0.134
Session            0.7783     0.6553   1.188
Drug1:Putamen_ki -86.0849   222.2221  -0.387

plot_model(MPH_Putamen_ConDivRatio, type = "pred", terms = c("Putamen_ki","Session","Drug"))

# Convergent/Total ratio score
MPH_Putamen_ConTotal <- lmer(Con_Total ~ 1 + Drug*Putamen_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_Putamen_ConTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Total ~ 1 + Drug * Putamen_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
   -54.1    -31.5     34.0    -68.1      179 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.0773 -0.4750  0.0937  0.5361  2.1925 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03480  0.1865  
 Residual             0.01858  0.1363  
Number of obs: 186, groups:  ID, 94

Fixed effects:
                  Estimate Std. Error t value
(Intercept)       0.510966   0.168119   3.039
Drug1             0.006385   0.077366   0.083
Putamen_ki        5.199394   9.757426   0.533
Session           0.032703   0.013777   2.374
Drug1:Putamen_ki -1.086851   4.556358  -0.239

plot_model(MPH_Putamen_ConTotal, type = "pred", terms = c("Putamen_ki","Session","Drug"))

# Convergent-Divergent difference score 
MPH_Putamen_ConDifference <- lmer(ConDiv_Dif ~ 1 + Drug*Putamen_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_Putamen_ConDifference), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: ConDiv_Dif ~ 1 + Drug * Putamen_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1609.9   1632.5   -797.9   1595.9      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.1311 -0.3807 -0.0473  0.4687  3.9547 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 232.2    15.24   
 Residual             135.0    11.62   
Number of obs: 188, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)         9.432     13.874   0.680
Drug1               1.538      6.593   0.233
Putamen_ki       -121.965    804.469  -0.152
Session             2.674      1.159   2.306
Drug1:Putamen_ki -123.085    388.327  -0.317

plot_model(MPH_Putamen_ConDifference, type = "pred", terms = c("Putamen_ki","Session","Drug"))

plot(MPH_Putamen_ConDifference)

qqnorm(residuals(MPH_Putamen_ConDifference))

# Divergent/Total ratio score
MPH_Putamen_DivTotal <- lmer(Div_Total ~ 1 + Drug*Putamen_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_Putamen_DivTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Div_Total ~ 1 + Drug * Putamen_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
   -53.9    -31.3     34.0    -67.9      179 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.19125 -0.53605 -0.09337  0.47557  2.07518 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03478  0.1865  
 Residual             0.01862  0.1364  
Number of obs: 186, groups:  ID, 94

Fixed effects:
                  Estimate Std. Error t value
(Intercept)       0.489172   0.168122   2.910
Drug1            -0.006116   0.077438  -0.079
Putamen_ki       -5.179289   9.757352  -0.531
Session          -0.032811   0.013789  -2.380
Drug1:Putamen_ki  1.074502   4.560597   0.236

plot_model(MPH_Putamen_DivTotal, type = "pred", terms = c("Putamen_ki","Session","Drug"))

## effects of SUL vs PBO below

# pure Convergent score
SUL_Putamen_PureConvergent <- lmer(Convergent_Pasta ~ 1 + Drug*Putamen_ki + Session + (1 | ID), data = df_SUL, REML=F) 
print(summary(SUL_Putamen_PureConvergent), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Convergent_Pasta ~ 1 + Drug * Putamen_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
  1495.8   1518.4   -740.9   1481.8      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.9680 -0.4589 -0.1037  0.3027  4.0126 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 178.33   13.354  
 Residual              57.99    7.615  
Number of obs: 188, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)       17.6047    11.3887   1.546
Drug1             -1.2096     4.2438  -0.285
Putamen_ki       -44.8624   669.5883  -0.067
Session            2.8012     0.7939   3.528
Drug1:Putamen_ki  53.4798   250.1742   0.214

plot_model(SUL_Putamen_PureConvergent, type = "pred", terms = c("Putamen_ki","Session","Drug"))

# Convergent/Divergent ratio score
SUL_Putamen_ConDivRatio <- lmer(Con_Div ~ 1 + Drug*Putamen_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_Putamen_ConDivRatio), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Div ~ 1 + Drug * Putamen_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
  1205.6   1227.9   -595.8   1191.6      171 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.8730 -0.3730 -0.2144  0.0592  5.7171 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 17.04    4.128   
 Residual             33.11    5.754   
Number of obs: 178, groups:  ID, 93

Fixed effects:
                  Estimate Std. Error t value
(Intercept)         7.1114     4.6993   1.513
Drug1               4.0429     3.3049   1.223
Putamen_ki       -232.5022   274.6479  -0.847
Session             0.7103     0.5885   1.207
Drug1:Putamen_ki -262.3971   194.9995  -1.346

plot_model(SUL_Putamen_ConDivRatio, type = "pred", terms = c("Putamen_ki","Session","Drug"))

# Convergent/Total ratio score
SUL_Putamen_ConTotal <- lmer(Con_Total ~ 1 + Drug*Putamen_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_Putamen_ConTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Total ~ 1 + Drug * Putamen_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
   -51.1    -28.4     32.5    -65.1      181 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.17868 -0.55058  0.07223  0.50396  1.82385 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03357  0.1832  
 Residual             0.01975  0.1405  
Number of obs: 188, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)       0.47846    0.16528   2.895
Drug1             0.01354    0.07831   0.173
Putamen_ki        2.82104    9.69880   0.291
Session           0.06784    0.01444   4.699
Drug1:Putamen_ki -1.88813    4.61658  -0.409

plot_model(SUL_Putamen_ConTotal, type = "pred", terms = c("Putamen_ki","Session","Drug"))

# Convergent-Divergent difference score 
SUL_Putamen_ConDifference <- lmer(ConDiv_Dif ~ 1 + Drug*Putamen_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_Putamen_ConDifference), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: ConDiv_Dif ~ 1 + Drug * Putamen_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
  1589.9   1612.6   -788.0   1575.9      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.7219 -0.4398 -0.0209  0.3943  3.8002 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 217.2    14.74   
 Residual             118.5    10.88   
Number of obs: 188, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)        6.5821    13.1780   0.499
Drug1             -0.7163     6.0655  -0.118
Putamen_ki       -15.8827   773.5282  -0.021
Session            3.2009     1.1208   2.856
Drug1:Putamen_ki   6.6629   357.5738   0.019

plot_model(SUL_Putamen_ConDifference, type = "pred", terms = c("Putamen_ki","Session","Drug"))

plot(SUL_Putamen_ConDifference)

qqnorm(residuals(SUL_Putamen_ConDifference))

# Divergent/Total ratio score
SUL_Putamen_DivTotal <- lmer(Div_Total ~ 1 + Drug*Putamen_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_Putamen_DivTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Div_Total ~ 1 + Drug * Putamen_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
   -51.1    -28.4     32.6    -65.1      181 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-1.82483 -0.50493 -0.07243  0.55023  2.17833 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03355  0.1832  
 Residual             0.01975  0.1405  
Number of obs: 188, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)       0.52141    0.16524   3.156
Drug1            -0.01341    0.07831  -0.171
Putamen_ki       -2.80731    9.69635  -0.290
Session          -0.06784    0.01444  -4.698
Drug1:Putamen_ki  1.87433    4.61667   0.406

plot_model(SUL_Putamen_DivTotal, type = "pred", terms = c("Putamen_ki","Session","Drug"))

############################# Test the effects of VS Ki #############################################################

## effects of MPH vs PBO below

# pure Convergent score
MPH_VS_PureConvergent <- lmer(Convergent_Pasta ~ 1 + Drug*VS_ki + Session + (1 | ID), data = df_MPH, REML=F) 
print(summary(MPH_VS_PureConvergent), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Convergent_Pasta ~ 1 + Drug * VS_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1520.0   1542.6   -753.0   1506.0      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.2164 -0.4128 -0.0839  0.3701  4.1958 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 163.41   12.783  
 Residual              77.03    8.777  
Number of obs: 188, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)   4.5639    11.5073   0.397
Drug1        -6.1717     5.0276  -1.228
VS_ki       907.3052   772.0248   1.175
Session       2.5662     0.8773   2.925
Drug1:VS_ki 435.6364   341.8961   1.274

plot_model(MPH_VS_PureConvergent, show.values=TRUE, value.offset = .3)

plot_model(MPH_VS_PureConvergent,  type = "pred", terms = c("VS_ki","Session","Drug"))

# Convergent/Divergent ratio score
MPH_VS_ConDivRatio <- lmer(Con_Div ~ 1 + Drug*VS_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_VS_ConDivRatio), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Div ~ 1 + Drug * VS_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1227.6   1249.7   -606.8   1213.6      168 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.7999 -0.3856 -0.2262 -0.0557  5.1621 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 20.07    4.480   
 Residual             43.21    6.574   
Number of obs: 175, groups:  ID, 91

Fixed effects:
            Estimate Std. Error t value
(Intercept)   1.7153     5.4730   0.313
Drug1        -1.5541     3.8669  -0.402
VS_ki       108.2781   359.4776   0.301
Session       0.8644     0.6529   1.324
Drug1:VS_ki 105.5954   262.7893   0.402

plot_model(MPH_VS_ConDivRatio, show.values=TRUE, value.offset = .3)

plot_model(MPH_VS_ConDivRatio, type = "pred", terms = c("VS_ki","Session","Drug"))

# Convergent/Total ratio score
MPH_VS_ConTotal <- lmer(Con_Total ~ 1 + Drug*VS_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_VS_ConTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Total ~ 1 + Drug * VS_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
   -55.1    -32.5     34.5    -69.1      179 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.1040 -0.5145  0.1039  0.5376  2.2060 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03467  0.1862  
 Residual             0.01846  0.1359  
Number of obs: 186, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)  0.47019    0.17019   2.763
Drug1       -0.08165    0.07822  -1.044
VS_ki        8.43062   11.39733   0.740
Session      0.03556    0.01371   2.593
Drug1:VS_ki  4.79022    5.32597   0.899

plot_model(MPH_VS_ConTotal, show.values=TRUE, value.offset = .3)

plot_model(MPH_VS_ConTotal, type = "pred", terms = c("VS_ki","Session","Drug"))

# Convergent-Divergent difference score 
MPH_VS_ConDifference <- lmer(ConDiv_Dif ~ 1 + Drug*VS_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_VS_ConDifference), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: ConDiv_Dif ~ 1 + Drug * VS_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
  1608.0   1630.7   -797.0   1594.0      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-3.1317 -0.3960 -0.0366  0.4012  3.9846 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 228.7    15.12   
 Residual             134.2    11.58   
Number of obs: 188, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)   -8.830     13.964  -0.632
Drug1         -5.837      6.634  -0.880
VS_ki       1078.450    934.423   1.154
Session        2.936      1.152   2.549
Drug1:VS_ki  363.473    451.169   0.806

plot_model(MPH_VS_ConDifference, type = "pred", terms = c("VS_ki","Session","Drug"))

plot(MPH_VS_ConDifference)

qqnorm(residuals(MPH_VS_ConDifference))

# Divergent/Total ratio score
MPH_VS_DivTotal <- lmer(Div_Total ~ 1 + Drug*VS_ki + Session + (1 | ID), data = df_MPH, REML=F)
print(summary(MPH_VS_DivTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Div_Total ~ 1 + Drug * VS_ki + Session + (1 | ID)
   Data: df_MPH

     AIC      BIC   logLik deviance df.resid 
   -54.9    -32.3     34.5    -68.9      179 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.2048 -0.5383 -0.1045  0.5148  2.1020 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03465  0.1861  
 Residual             0.01850  0.1360  
Number of obs: 186, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)  0.53043    0.17019   3.117
Drug1        0.08178    0.07829   1.045
VS_ki       -8.44068   11.39708  -0.741
Session     -0.03566    0.01373  -2.598
Drug1:VS_ki -4.79473    5.33086  -0.899

plot_model(MPH_VS_DivTotal, show.values=TRUE, value.offset = .3)

plot_model(MPH_VS_DivTotal, type = "pred", terms = c("VS_ki","Session","Drug"))

## effects of SUL vs PBO below

# pure Convergent score
SUL_VS_PureConvergent <- lmer(Convergent_Pasta ~ 1 + Drug*VS_ki + Session + (1 | ID), data = df_SUL, REML=F) 
print(summary(SUL_VS_PureConvergent), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Convergent_Pasta ~ 1 + Drug * VS_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
  1493.9   1516.6   -740.0   1479.9      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.9575 -0.4389 -0.0915  0.3183  4.0199 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 176.62   13.290  
 Residual              57.41    7.577  
Number of obs: 188, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)   5.9742    11.4791   0.520
Drug1        -4.5189     4.2738  -1.057
VS_ki       752.8912   779.0124   0.966
Session       2.7578     0.7898   3.492
Drug1:VS_ki 289.3428   290.9810   0.994

plot_model(SUL_VS_PureConvergent, show.values=TRUE, value.offset = .3)

plot_model(SUL_VS_PureConvergent, type = "pred", terms = c("VS_ki","Session","Drug"))

# Convergent/Divergent ratio score
SUL_VS_ConDivRatio <- lmer(Con_Div ~ 1 + Drug*VS_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_VS_ConDivRatio), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Div ~ 1 + Drug * VS_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
  1207.0   1229.3   -596.5   1193.0      171 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.8475 -0.3906 -0.1878  0.0142  5.7288 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 16.72    4.089   
 Residual             33.68    5.803   
Number of obs: 178, groups:  ID, 93

Fixed effects:
             Estimate Std. Error t value
(Intercept)    6.0577     4.8061   1.260
Drug1          2.4641     3.3765   0.730
VS_ki       -193.9864   322.5798  -0.601
Session        0.6976     0.5918   1.179
Drug1:VS_ki -193.7748   229.5201  -0.844

plot_model(SUL_VS_ConDivRatio, show.values=TRUE, value.offset = .3)

plot_model(SUL_VS_ConDivRatio, type = "pred", terms = c("VS_ki","Session","Drug"))

# Convergent/Total ratio score
SUL_VS_ConTotal <- lmer(Con_Total ~ 1 + Drug*VS_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_VS_ConTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Con_Total ~ 1 + Drug * VS_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
   -51.3    -28.7     32.7    -65.3      181 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.15283 -0.54150  0.06619  0.52389  1.83639 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03365  0.1834  
 Residual             0.01968  0.1403  
Number of obs: 188, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)  0.54692    0.16747   3.266
Drug1        0.03660    0.07912   0.463
VS_ki       -1.51742   11.34328  -0.134
Session      0.06837    0.01441   4.744
Drug1:VS_ki -3.76739    5.38705  -0.699

plot_model(SUL_VS_ConTotal, show.values=TRUE, value.offset = .3)

plot_model(SUL_VS_ConTotal, type = "pred", terms = c("VS_ki","Session","Drug"))

# Convergent-Divergent difference score 
SUL_VS_ConDifference <- lmer(ConDiv_Dif ~ 1 + Drug*VS_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_VS_ConDifference), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: ConDiv_Dif ~ 1 + Drug * VS_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
  1589.1   1611.8   -787.6   1575.1      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-2.7051 -0.4561 -0.0262  0.3903  3.8065 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 214.9    14.66   
 Residual             118.4    10.88   
Number of obs: 188, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)   -5.174     13.292  -0.389
Drug1         -1.919      6.138  -0.313
VS_ki        796.974    900.504   0.885
Session        3.143      1.120   2.807
Drug1:VS_ki   90.778    417.924   0.217

plot_model(SUL_VS_ConDifference, type = "pred", terms = c("VS_ki","Session","Drug"))

plot(SUL_VS_ConDifference)

qqnorm(residuals(SUL_VS_ConDifference))

# Divergent/Total ratio score
SUL_VS_DivTotal <- lmer(Div_Total ~ 1 + Drug*VS_ki + Session + (1 | ID), data = df_SUL, REML=F)
print(summary(SUL_VS_DivTotal), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Div_Total ~ 1 + Drug * VS_ki + Session + (1 | ID)
   Data: df_SUL

     AIC      BIC   logLik deviance df.resid 
   -51.4    -28.7     32.7    -65.4      181 

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-1.83744 -0.52377 -0.06621  0.54220  2.15188 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.03362  0.1834  
 Residual             0.01968  0.1403  
Number of obs: 188, groups:  ID, 94

Fixed effects:
            Estimate Std. Error t value
(Intercept)  0.45326    0.16743   2.707
Drug1       -0.03678    0.07912  -0.465
VS_ki        1.51227   11.34036   0.133
Session     -0.06836    0.01441  -4.744
Drug1:VS_ki  3.77236    5.38706   0.700

plot_model(SUL_VS_DivTotal, show.values=TRUE, value.offset = .3)

plot_model(SUL_VS_DivTotal, type = "pred", terms = c("VS_ki","Session","Drug"))

NA
NA

# Factorize the data into score types and test the effects of score type as a factor on scores (1=convergent, 2=divergent)

FactorizedScores <- read.csv('~/Dropbox/DCCN/Creativity/CreativityIDsheet_MissingDataRemoved_ScoresFactorized.csv')

# turn drug conditions into factor levels
FactorizedScores$Drug <- factor(FactorizedScores$Drug, levels = c("MPH","SUL","PBO"));
FactorizedScores$ScoreType <- factor(FactorizedScores$ScoreType);
hs <- function(x) {
    y <- 0.5*exp(-x)*(exp(2*x)-1)
    return(y)
}  # inverse sine transormation
#FactorizedScores$Score <- hs(FactorizedScores$Score); # zscored the scores
FactorizedScores$Score <-scale(FactorizedScores$Score, center = TRUE, scale = TRUE)

# ScoreType 1 refers to convergent, 2 refers to divergent

# set contrasts to sum-to-zero
options(contrasts=c("contr.sum", "contr.poly"))

# set two dataframes for the contrast between MPH and PBO - between SUL and PBO
df_MPH_Factorized <- dplyr::filter(FactorizedScores, Drug %in% c("MPH","PBO"))  
df_SUL_Factorized <- dplyr::filter(FactorizedScores, Drug %in% c("SUL","PBO"))

## effects of MPH vs PBO on factorized model

MPH_VS_FactorizedScore <- lmer(Score ~ 1 + Drug*VS_ki*ScoreType + Session + (1 | ID), data = df_MPH_Factorized, REML=F) 
print(summary(MPH_VS_FactorizedScore), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Score ~ 1 + Drug * VS_ki * ScoreType + Session + (1 | ID)
   Data: df_MPH_Factorized

     AIC      BIC   logLik deviance df.resid 
   996.0   1039.2   -487.0    974.0      365 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.9628 -0.6380 -0.1269  0.4391  5.9261 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.09562  0.3092  
 Residual             0.70010  0.8367  
Number of obs: 376, groups:  ID, 94

Fixed effects:
                       Estimate Std. Error t value
(Intercept)            -0.46831    0.43371  -1.080
Drug1                  -0.20533    0.33832  -0.607
VS_ki                  25.82118   28.28359   0.913
ScoreType1             -0.07783    0.33367  -0.233
Session                 0.05047    0.05560   0.908
Drug1:VS_ki            16.39899   23.01029   0.713
Drug1:ScoreType1       -0.10384    0.33367  -0.311
VS_ki:ScoreType1       36.36198   22.72164   1.600
Drug1:VS_ki:ScoreType1  6.17207   22.72164   0.272

plot_model(MPH_VS_FactorizedScore, show.values=TRUE, value.offset = .3)

plot_model(MPH_VS_FactorizedScore,  type = "pred", terms = c("VS_ki","ScoreType","Drug"))

plot(MPH_VS_FactorizedScore)

qqnorm(residuals(MPH_VS_FactorizedScore))

MPH_Putamen_FactorizedScore <- lmer(Score ~ 1 + Drug*Putamen_ki*ScoreType + Session + (1 | ID), data = df_MPH_Factorized, REML=F) 
print(summary(MPH_Putamen_FactorizedScore), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Score ~ 1 + Drug * Putamen_ki * ScoreType + Session + (1 | ID)
   Data: df_MPH_Factorized

     AIC      BIC   logLik deviance df.resid 
   999.5   1042.7   -488.7    977.5      365 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.8945 -0.6345 -0.1318  0.4171  6.0511 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.09574  0.3094  
 Residual             0.70728  0.8410  
Number of obs: 376, groups:  ID, 94

Fixed effects:
                             Estimate Std. Error t value
(Intercept)                  -0.05924    0.42951  -0.138
Drug1                        -0.01969    0.33673  -0.058
Putamen_ki                   -1.14636   24.26996  -0.047
ScoreType1                    0.55270    0.33113   1.669
Session                       0.04335    0.05601   0.774
Drug1:Putamen_ki              3.18561   19.83745   0.161
Drug1:ScoreType1              0.16026    0.33113   0.484
Putamen_ki:ScoreType1        -6.01111   19.53336  -0.308
Drug1:Putamen_ki:ScoreType1 -10.36694   19.53336  -0.531

plot_model(MPH_Putamen_FactorizedScore, show.values=TRUE, value.offset = .3)

plot_model(MPH_Putamen_FactorizedScore,  type = "pred", terms = c("Putamen_ki","ScoreType","Drug"))

plot(MPH_Putamen_FactorizedScore)

qqnorm(residuals(MPH_Putamen_FactorizedScore))

MPH_Caudate_FactorizedScore <- lmer(Score ~ 1 + Drug*Caudate_ki*ScoreType + Session + (1 | ID), data = df_MPH_Factorized, REML=F) 
print(summary(MPH_Caudate_FactorizedScore), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Score ~ 1 + Drug * Caudate_ki * ScoreType + Session + (1 | ID)
   Data: df_MPH_Factorized

     AIC      BIC   logLik deviance df.resid 
   999.4   1042.6   -488.7    977.4      365 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.9110 -0.6321 -0.1090  0.4160  6.0466 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.09567  0.3093  
 Residual             0.70708  0.8409  
Number of obs: 376, groups:  ID, 94

Fixed effects:
                             Estimate Std. Error t value
(Intercept)                  0.008465   0.392731   0.022
Drug1                       -0.160904   0.301161  -0.534
Caudate_ki                  -6.303045  26.389338  -0.239
ScoreType1                   0.454306   0.300111   1.514
Session                      0.043955   0.055316   0.795
Drug1:Caudate_ki            13.923455  21.293132   0.654
Drug1:ScoreType1            -0.074264   0.300111  -0.247
Caudate_ki:ScoreType1       -0.188328  21.233979  -0.009
Drug1:Caudate_ki:ScoreType1  4.311609  21.233979   0.203

plot_model(MPH_Caudate_FactorizedScore, show.values=TRUE, value.offset = .3)

plot_model(MPH_Caudate_FactorizedScore,  type = "pred", terms = c("Caudate_ki","ScoreType","Drug"))

plot(MPH_Caudate_FactorizedScore)

qqnorm(residuals(MPH_Caudate_FactorizedScore))

## effects of SUL vs PBO on factorized model


SUL_VS_FactorizedScore <- lmer(Score ~ 1 + Drug*VS_ki*ScoreType + Session + (1 | ID), data = df_SUL_Factorized, REML=F) 
print(summary(SUL_VS_FactorizedScore), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Score ~ 1 + Drug * VS_ki * ScoreType + Session + (1 | ID)
   Data: df_SUL_Factorized

     AIC      BIC   logLik deviance df.resid 
   972.1   1015.3   -475.1    950.1      365 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.7721 -0.6374 -0.1463  0.4634  4.2719 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.1462   0.3824  
 Residual             0.6207   0.7878  
Number of obs: 376, groups:  ID, 94

Fixed effects:
                       Estimate Std. Error t value
(Intercept)            -0.58791    0.44416  -1.324
Drug1                  -0.25778    0.31422  -0.820
VS_ki                  24.68895   29.91420   0.825
ScoreType1             -0.03210    0.31418  -0.102
Session                 0.10353    0.05470   1.893
Drug1:VS_ki            17.51553   21.39430   0.819
Drug1:ScoreType1       -0.05811    0.31418  -0.185
VS_ki:ScoreType1       33.65525   21.39408   1.573
Drug1:VS_ki:ScoreType1  3.46534   21.39408   0.162

plot_model(SUL_VS_FactorizedScore, show.values=TRUE, value.offset = .3)

plot_model(SUL_VS_FactorizedScore,  type = "pred", terms = c("VS_ki","ScoreType","Drug"))

plot(SUL_VS_FactorizedScore)

qqnorm(residuals(SUL_VS_FactorizedScore))

SUL_Putamen_FactorizedScore <- lmer(Score ~ 1 + Drug*Putamen_ki*ScoreType + Session + (1 | ID), data = df_SUL_Factorized, REML=F) 
print(summary(SUL_Putamen_FactorizedScore), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Score ~ 1 + Drug * Putamen_ki * ScoreType + Session + (1 | ID)
   Data: df_SUL_Factorized

     AIC      BIC   logLik deviance df.resid 
   975.8   1019.1   -476.9    953.8      365 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.8211 -0.6335 -0.1333  0.4659  4.2397 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.1463   0.3825  
 Residual             0.6276   0.7922  
Number of obs: 376, groups:  ID, 94

Fixed effects:
                              Estimate Std. Error t value
(Intercept)                 -0.1788738  0.4399031  -0.407
Drug1                       -0.0659740  0.3121581  -0.211
Putamen_ki                  -3.5193419 25.6718981  -0.137
ScoreType1                   0.3914818  0.3119329   1.255
Session                      0.1083277  0.0550036   1.969
Drug1:Putamen_ki             3.7289515 18.4033478   0.203
Drug1:ScoreType1            -0.0009515  0.3119329  -0.003
Putamen_ki:ScoreType1        3.9571244 18.4008205   0.215
Drug1:Putamen_ki:ScoreType1 -0.3987063 18.4008205  -0.022

plot_model(SUL_Putamen_FactorizedScore, show.values=TRUE, value.offset = .3)

plot_model(SUL_Putamen_FactorizedScore,  type = "pred", terms = c("Putamen_ki","ScoreType","Drug"))

plot(SUL_Putamen_FactorizedScore)

qqnorm(residuals(SUL_Putamen_FactorizedScore))

SUL_Caudate_FactorizedScore <- lmer(Score ~ 1 + Drug*Caudate_ki*ScoreType + Session + (1 | ID), data = df_SUL_Factorized, REML=F) 
print(summary(SUL_Caudate_FactorizedScore), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Score ~ 1 + Drug * Caudate_ki * ScoreType + Session + (1 | ID)
   Data: df_SUL_Factorized

     AIC      BIC   logLik deviance df.resid 
   975.2   1018.4   -476.6    953.2      365 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.9221 -0.6350 -0.1389  0.4458  4.2110 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.1464   0.3826  
 Residual             0.6264   0.7914  
Number of obs: 376, groups:  ID, 94

Fixed effects:
                            Estimate Std. Error t value
(Intercept)                 -0.09991    0.40323  -0.248
Drug1                       -0.15658    0.28257  -0.554
Caudate_ki                  -9.66734   27.83538  -0.347
ScoreType1                   0.37966    0.28246   1.344
Session                      0.10688    0.05489   1.947
Drug1:Caudate_ki            10.97217   20.01212   0.548
Drug1:ScoreType1            -0.14891    0.28246  -0.527
Caudate_ki:ScoreType1        5.60054   19.98541   0.280
Drug1:Caudate_ki:ScoreType1 10.10048   19.98541   0.505

plot_model(SUL_Caudate_FactorizedScore, show.values=TRUE, value.offset = .3)

plot_model(SUL_Caudate_FactorizedScore,  type = "pred", terms = c("Caudate_ki","ScoreType","Drug"))

plot(SUL_Caudate_FactorizedScore)

qqnorm(residuals(SUL_Caudate_FactorizedScore))

## effects of Ki during placebo on factorized model

OnlyPBO <- dplyr::filter(FactorizedScores, Drug %in% c("PBO"))  

PBO_VS_FactorizedScore <- lmer(Score ~ 1 + VS_ki*ScoreType + Session + (1 | ID), data = OnlyPBO, REML=F) 

boundary (singular) fit: see ?isSingular

print(summary(PBO_VS_FactorizedScore), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Score ~ 1 + VS_ki * ScoreType + Session + (1 | ID)
   Data: OnlyPBO

     AIC      BIC   logLik deviance df.resid 
   503.5    526.2   -244.8    489.5      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.7456 -0.6122 -0.1823  0.4128  4.3174 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.0000   0.0000  
 Residual             0.7913   0.8896  
Number of obs: 188, groups:  ID, 94

Fixed effects:
                 Estimate Std. Error t value
(Intercept)      -0.25688    0.51192  -0.502
VS_ki             9.62638   34.33228   0.280
ScoreType1        0.02601    0.50169   0.052
Session           0.04565    0.08047   0.567
VS_ki:ScoreType1 30.18991   34.16245   0.884
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see ?isSingular

plot_model(PBO_VS_FactorizedScore, show.values=TRUE, value.offset = .3)

plot_model(PBO_VS_FactorizedScore,  type = "pred", terms = c("VS_ki","ScoreType"))

plot(PBO_VS_FactorizedScore)

qqnorm(residuals(PBO_VS_FactorizedScore))

PBO_Putamen_FactorizedScore <- lmer(Score ~ 1 + Putamen_ki*ScoreType + Session + (1 | ID), data = OnlyPBO, REML=F) 

boundary (singular) fit: see ?isSingular

print(summary(PBO_Putamen_FactorizedScore), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Score ~ 1 + Putamen_ki * ScoreType + Session + (1 | ID)
   Data: OnlyPBO

     AIC      BIC   logLik deviance df.resid 
   504.3    527.0   -245.2    490.3      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.6556 -0.6098 -0.1839  0.4319  4.2902 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.0000   0.0000  
 Residual             0.7947   0.8915  
Number of obs: 188, groups:  ID, 94

Fixed effects:
                      Estimate Std. Error t value
(Intercept)           -0.04644    0.50472  -0.092
Putamen_ki            -4.60570   29.50633  -0.156
ScoreType1             0.39243    0.49640   0.791
Session                0.04945    0.08086   0.612
Putamen_ki:ScoreType1  4.35583   29.28233   0.149
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see ?isSingular

plot_model(PBO_Putamen_FactorizedScore, show.values=TRUE, value.offset = .3)

plot_model(PBO_Putamen_FactorizedScore,  type = "pred", terms = c("Putamen_ki","ScoreType"))

plot(PBO_Putamen_FactorizedScore)

qqnorm(residuals(PBO_Putamen_FactorizedScore))

PBO_Caudate_FactorizedScore <- lmer(Score ~ 1 + Caudate_ki*ScoreType + Session + (1 | ID), data = OnlyPBO, REML=F) 

boundary (singular) fit: see ?isSingular

print(summary(PBO_Caudate_FactorizedScore), corr=F) # print summary without fixed effect correlation matrix

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: Score ~ 1 + Caudate_ki * ScoreType + Session + (1 | ID)
   Data: OnlyPBO

     AIC      BIC   logLik deviance df.resid 
   503.9    526.6   -245.0    489.9      181 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-1.6763 -0.6183 -0.2027  0.3962  4.2602 

Random effects:
 Groups   Name        Variance Std.Dev.
 ID       (Intercept) 0.0000   0.0000  
 Residual             0.7931   0.8906  
Number of obs: 188, groups:  ID, 94

Fixed effects:
                       Estimate Std. Error t value
(Intercept)             0.16080    0.47189   0.341
Caudate_ki            -20.25790   31.80862  -0.637
ScoreType1              0.52857    0.44951   1.176
Session                 0.04874    0.08017   0.608
Caudate_ki:ScoreType1  -4.49994   31.80427  -0.141
optimizer (nloptwrap) convergence code: 0 (OK)
boundary (singular) fit: see ?isSingular

plot_model(PBO_Caudate_FactorizedScore, show.values=TRUE, value.offset = .3)

plot_model(PBO_Caudate_FactorizedScore,  type = "pred", terms = c("Caudate_ki","ScoreType"))

plot(PBO_Caudate_FactorizedScore)

qqnorm(residuals(PBO_Caudate_FactorizedScore))

NA
NA

Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.